##################################################
## Authors: Stephanie L. DeMora, Jennifer L. Merolla, Maricruz Ariana Osorio, Christian Lindke, and Sean Long
## Project: Emotional Engagement in the Year of the Woman
## Purpose: Replication File for Publication: In-Text Tables & Figures
## Date: October 27, 2025
##################################################

# Load in all packages
Packages <- c("here", "car", "survey", "margins", "interplot","aod","dplyr","foreign","nnet","ggplot2","reshape2","Rmisc",
              "tidyverse", "broom","dotwhisker","janitor","qualtRics","jtools","mediation","gtools","gridExtra", "huxtable", 
              "DT", "ggsci", "hrbrthemes", "haven", "effectsize", "readr", "purrr", "emmeans")


lapply(Packages, library, character.only = TRUE)
options(scipen=999)

# Set plot theme:
steph_theme <- list(
  hrbrthemes::theme_ipsum_rc(base_size = 12,
                             grid= "",
                             plot_title_size = 20,axis_text_size = 18),
  ggsci::scale_color_npg(),
  ggsci::scale_fill_npg(),
  theme(legend.position = "top",
        plot.caption = element_text(hjust = 0, face= "italic"), #Default is hjust=1
        plot.title.position = "plot", #NEW parameter. Apply for subtitle too.
        plot.caption.position =  "plot"))

# GRAYSCALE Plot Theme:
steph_theme_gray <- list(
  hrbrthemes::theme_ipsum_rc(base_size = 12,
                             grid= "",
                             plot_title_size = 20,axis_text_size = 18),
  ggsci::scale_color_material("grey"),
  ggsci::scale_fill_material("grey"),
  theme(legend.position = "top",
        plot.caption = element_text(hjust = 0, face= "italic"), #Default is hjust=1
        plot.title.position = "plot", #NEW parameter. Apply for subtitle too.
        plot.caption.position =  "plot"))

# Read in data:
df <- readRDS("Qualtrics_2019.RDS")
df %>%
  filter(TREATMENT == "Trump" | TREATMENT == "Clinton" | TREATMENT == "Alexa") -> df
df$TREATMENT_sort <- factor(df$TREATMENT, levels = c("Alexa","Trump","Clinton"))

df_TESS <- readRDS("NORC.RDS")
df_TESS$TREATMENT_sort <- factor(df_TESS$TREATMENT, levels = c("Alexa","Trump","Clinton"))
df_TESS %>%
  filter(TREATMENT == "Trump" | TREATMENT == "Clinton" | TREATMENT == "Alexa") -> df_TESS
df_TESS$TREATMENT_sort <- factor(df_TESS$TREATMENT, levels = c("Alexa","Trump","Clinton"))

anes <- readRDS("ANES16.RDS")
anes18 <- readRDS("ANES18.RDS")

# ---------------------------------------------------------------------------------------------------------------
# Figure 1: Mean Level of Anger by Experimental Condition, Study 1
# ---------------------------------------------------------------------------------------------------------------
options(scipen=999)
df1 <- subset(df,df$anger != "<NA>")
summary(lm(anger ~ TREATMENT, data = df))
summarySE(df1, measurevar="anger", groupvars=c("TREATMENT"), conf.interval = .90)
TukeyHSD(aov(anger~TREATMENT,data = df))

tgc <- summarySE(df1, measurevar="anger", groupvars=c("TREATMENT"), conf.interval = .90)
tgc$TREATMENT <- c("Control","Clinton","Trump")
tgc$low = tgc$anger-tgc$ci; tgc$high = tgc$anger+tgc$ci; tgc$sort <- c(3, 2, 1)
tgc %>%
  arrange(sort) -> tgc

ggplot(data=tgc, aes(x=reorder(TREATMENT, sort), y=anger)) +
  geom_bar(position = "dodge", stat="identity") +
  geom_errorbar(position = "dodge",aes(ymin=low, ymax=high)) + 
  steph_theme_gray + 
  labs(title="Mean Levels of Anger", 
       x="Treatment Condition", y = "Average level of anger", caption = "Fall 2019 Qualtrics Survey",subtitle = "By Experimental Condition", fill = "Gender:") + 
  theme(legend.position = "none")

# ---------------------------------------------------------------------------------------------------------------
# In Text Only:

# Replication of Table 1 but controlling for pre-treatment vote choice on ANGER:
anger1 <- lm(anger ~ TREATMENT, data = df)
anger2 <- lm(anger ~ TREATMENT + Vote_Choice, data = df)
anger3 <- lm(anger ~ TREATMENT + Vote_Choice + party_leanersinc + GENDER + EDUCATION_r, data = df)
export_summs(anger1,anger2,anger3, error_format = "({p.value})")

# Replication of Table 1 but controlling for pre-treatment vote choice on ENTHUSIASM:
anger1 <- lm(hope ~ TREATMENT, data = df)
anger2 <- lm(hope ~ TREATMENT + Vote_Choice, data = df)
anger3 <- lm(hope ~ TREATMENT + Vote_Choice + party_leanersinc + GENDER + EDUCATION_r, data = df)
export_summs(anger1,anger2,anger3, error_format = "({p.value})")

# Replication of Table 1 but controlling for pre-treatment vote choice on ANXIETY:
anger1 <- lm(fear ~ TREATMENT, data = df)
anger2 <- lm(fear ~ TREATMENT + Vote_Choice, data = df)
anger3 <- lm(fear ~ TREATMENT + Vote_Choice + party_leanersinc + GENDER + EDUCATION_r, data = df)
export_summs(anger1,anger2,anger3, error_format = "({p.value})")


# ---------------------------------------------------------------------------------------------------------------
# Figure 2. Change in Anger Moving from the Control Group to the Trump Condition, by Gender and Party, Qualtrics Panel, Study 1
# ---------------------------------------------------------------------------------------------------------------
# calculate margins at 90% instead of 95%. Cures issue with original function.
margins_summary <- function(model, ..., level = 0.90, by_factor = TRUE) {
  mar <- margins(model, ...)
  summary(mar, level = 0.90, by_factor = by_factor)
}

# Republicans & Democrats
df1 <- subset(df,df$party_leanersinc!='NA')
anger1 <- lm(anger ~ TREATMENT*party_leanersinc, data = df1)
margins_summary(anger1, at = list(party_leanersinc = range(df1$party_leanersinc)))[c(3:6),]

# MEN AND WOMEN
df1<-subset(df,df$GENDER!='NA')
anger1 <- lm(anger ~ TREATMENT*GENDER, data = df1)
margins_summary(anger1, at = list(GENDER = range(df1$GENDER)))[c(3:6),]

# Republican & Democratic Women:
df1<-subset(df,df$GENDER != 'NA')
df1<-subset(df1,df1$party_leanersinc!='NA')
anger1 <- lm(anger ~ TREATMENT*party_leanersinc*GENDER, data = df1)
data1 <- margins_summary(anger1, at = list(party_leanersinc = range(df1$party_leanersinc), GENDER = range(df1$GENDER)))[c(13:16),]

data1 %>%
  dplyr::mutate(GENDER = recode(GENDER,
                                `0` = "Men",
                                `1` = "Women")) %>%
  dplyr::mutate(party_leanersinc = recode(party_leanersinc,
                                          `0` = "Republican",
                                          `1` = "Democratic")) %>%
  dplyr::mutate(term = paste(party_leanersinc, GENDER, sep = " ")) %>%
  dplyr::rename(estimate = AME,
                conf.low = lower,
                conf.high = upper) %>%
  dplyr::select(-GENDER, -party_leanersinc, -SE, -z, -p, -factor) -> data1

# Men and Women:
df1<-subset(df,df$GENDER != "NA")
anger1 <- lm(anger ~ TREATMENT*GENDER, data = df1)
data <- margins_summary(anger1, at = list(GENDER = range(df1$GENDER)))[5:6,]

data %>%
  dplyr::mutate(term = recode(GENDER,
                              `0` = "Men",
                              `1` = "Women")) %>%
  dplyr::rename(estimate = AME,
                conf.low = lower,
                conf.high = upper) %>%
  dplyr::select(-SE, -z, -p, -factor, -GENDER) -> data

data <- rbind(data,data1)

dwplot(data,
       dot_args = list(color = "black"), # color for the dot
       whisker_args = list(color = "grey")) +
  theme_bw() +
  steph_theme +
  theme(legend.title = element_blank(), legend.position="none") + 
  xlab("") +
  geom_vline(xintercept = 0, colour = "grey60", linetype = 2) -> Fig2; Fig2

# ---------------------------------------------------------------------------------------------------------------
# In Text Analyses: Education

df1<-subset(df,df$GENDER != 'NA')
df1<-subset(df1,df1$EDUCATION_r!='NA')
anger1 <- lm(anger ~ TREATMENT*GENDER*EDUCATION_r, data = df1)
data2 <- margins_summary(anger1, at = list(GENDER = range(df1$GENDER), EDUCATION_r = range(df1$EDUCATION_r)))[c(9:16),]

data2 %>%
  dplyr::mutate(GENDER = recode(GENDER,
                                `0` = "Men",
                                `1` = "Women")) %>%
  dplyr::mutate(EDUCATION_r = recode(EDUCATION_r,
                                     `0` = "Less Educated",
                                     `1` = "Highly Educated")) %>%
  dplyr::mutate(Treatment = recode(factor,
                                   `TREATMENTClinton` = "Clinton",
                                   `TREATMENTTrump` = "Trump")) %>%
  #filter(GENDER == "Women" & party_leanersinc == "Democratic") %>%
  dplyr::mutate(term = paste(EDUCATION_r, GENDER, sep = " ")) %>%
  dplyr::rename(estimate = AME,
                conf.low = lower,
                conf.high = upper) %>%
  dplyr::select(-GENDER, -EDUCATION_r, -SE, -z, -factor) -> data2; data2


# ---------------------------------------------------------------------------------------------------------------
# Figure 3: Mean Level of Anger by Experimental Condition, NORC Panel, Study 2
# ---------------------------------------------------------------------------------------------------------------

df_TESS1 <- subset(df_TESS,df_TESS$anger != 'NA')
tgc <- summarySE(df_TESS1, measurevar="anger", groupvars=c("TREATMENT_sort"), conf.interval = .90)
tgc$low = tgc$anger-tgc$ci; tgc$high = tgc$anger+tgc$ci; tgc$sort <- c(3,1,2)
tgc %>%
  arrange(sort) -> tgc
tgc$TREATMENT_sort <- c("Trump","Clinton","Control")

ggplot(data=tgc, aes(x=reorder(TREATMENT_sort, sort), y=anger)) +
  geom_bar(position = "dodge", stat="identity") +
  geom_errorbar(position = "dodge",aes(ymin=low, ymax=high)) + 
  steph_theme + 
  labs(title="Mean Levels of Anger", 
       x="Treatment Condition", y = "Average level of anger", caption = "Fall 2019 NORC Survey",subtitle = "By Experimental Condition", fill = "Gender:") + 
  theme(legend.position = "none")

# ---------------------------------------------------------------------------------------------------------------
# Figure 4: Change in Anger Moving from the Control Group to the Trump and Clinton Conditions, by Party, NORC, Study 2
# ---------------------------------------------------------------------------------------------------------------
df_TESS %>%
  filter(complete.cases(party, TREATMENT, anger)) -> df_TESS1

df_TESS1$TREATMENT <- droplevels(df_TESS1$TREATMENT)
anger1 <- lm(anger ~ TREATMENT*party, data = df_TESS1)
data <- margins_summary(anger1, at = list(party = range(df_TESS1$party)))[3:6,]
data %>%
  dplyr::mutate(model = recode(party,
                               `0` = "Republican Women",
                               `1` = "Democratic Women")) %>%
  dplyr::rename(estimate = AME,
                conf.low = lower,
                conf.high = upper,
                term = factor) %>%
  dplyr::select(-party, -SE, -z, -p) %>%
  mutate(term = str_sub(term, 10)) -> data

dwplot(data) +
  theme_bw() +
  steph_theme +
  xlab("") +
  theme(legend.justification=c(.02, .993), legend.position=c(.02, .99),
        legend.title = element_blank(), legend.background =
          element_rect(color="gray90")) +
  scale_colour_grey(start = .3, end = .7,
                    name = "model") +
  geom_vline(xintercept = 0, colour = "grey60", linetype = 2)


# ---------------------------------------------------------------------------------------------------------------
# Figure 5: Mean Enthusiasm by Experimental Condition 
# ---------------------------------------------------------------------------------------------------------------
df_TESS1 <- subset(df_TESS,df_TESS$hope != 'NA')
tgc <- summarySE(df_TESS1, measurevar="hope", groupvars=c("TREATMENT"), conf.interval = .90)
tgc$low = tgc$hope-tgc$ci; tgc$high = tgc$hope+tgc$ci; tgc$sort <- c(3, 2, 1)
tgc %>%
  arrange(sort) -> tgc
tgc$TREATMENT <- c("Trump","Clinton","Control")

enthus2_plot <- ggplot(data=tgc, aes(x=reorder(TREATMENT, sort), y=hope)) +
  geom_bar(position = "dodge", stat="identity") +
  geom_errorbar(position = "dodge",aes(ymin=low, ymax=high)) + 
  steph_theme_gray + 
  labs(title="", 
       x="Treatment Condition", y = "Average level of enthusiasm", caption = "Fall 2019 NORC Survey",subtitle = "By Experimental Condition", fill = "Gender:") + 
  theme(legend.position = "none")

df1 <- subset(df,df$hope != "<NA>")
tgc <- summarySE(df1, measurevar="hope", groupvars=c("TREATMENT"), conf.interval = .90)
tgc$low = tgc$hope-tgc$ci; tgc$high = tgc$hope+tgc$ci; tgc$sort <- c(3, 2, 1)
tgc %>%
  arrange(sort) -> tgc
tgc$TREATMENT <- c("Trump","Clinton","Control")

enthus1_plot <- ggplot(data=tgc, aes(x=reorder(TREATMENT, sort), y=hope)) +
  geom_bar(position = "dodge", stat="identity") +
  geom_errorbar(position = "dodge",aes(ymin=low, ymax=high)) + 
  steph_theme_gray + 
  labs(title="Mean Levels of Enthusiasm", 
       x="Treatment Condition", y = "Average level of enthusiasm", caption = "Fall 2019 Qualtrics Survey",subtitle = "By Experimental Condition", fill = "Gender:") + 
  theme(legend.position = "none")

gridExtra::grid.arrange(enthus1_plot, enthus2_plot, ncol = 2)



# ---------------------------------------------------------------------------------------------------------------
# In Text Analyses: Education

df1<-subset(df,df$GENDER != 'NA')
df1<-subset(df1,df1$EDUCATION_r!='NA')
anger1 <- lm(hope ~ TREATMENT*GENDER*EDUCATION_r, data = df1)
data2 <- margins_summary(anger1, at = list(GENDER = range(df1$GENDER), EDUCATION_r = range(df1$EDUCATION_r)))[c(9:16),]

data2 %>%
  dplyr::mutate(GENDER = recode(GENDER,
                                `0` = "Men",
                                `1` = "Women")) %>%
  dplyr::mutate(EDUCATION_r = recode(EDUCATION_r,
                                     `0` = "Less Educated",
                                     `1` = "Highly Educated")) %>%
  dplyr::mutate(Treatment = recode(factor,
                                   `TREATMENTClinton` = "Clinton",
                                   `TREATMENTTrump` = "Trump")) %>%
  #filter(GENDER == "Women" & party_leanersinc == "Democratic") %>%
  dplyr::mutate(term = paste(EDUCATION_r, GENDER, sep = " ")) %>%
  dplyr::rename(estimate = AME,
                conf.low = lower,
                conf.high = upper) %>%
  dplyr::select(-GENDER, -EDUCATION_r, -SE, -z, -factor) -> data2; data2


# ---------------------------------------------------------------------------------------------------------------
# Figure 6: Change in Anger Moving from the Control Group to the Trump and Clinton Conditions, by Party, NORC, Study 2
# ---------------------------------------------------------------------------------------------------------------

df_TESS %>%
  filter(complete.cases(party, TREATMENT, hope)) -> df_TESS1

df_TESS1$TREATMENT <- droplevels(df_TESS1$TREATMENT)
anger1 <- lm(hope ~ TREATMENT*party, data = df_TESS1)
data <- margins_summary(anger1, at = list(party = range(df_TESS1$party)))[3:6,]
data %>%
  dplyr::mutate(model = recode(party,
                               `0` = "Republican Women",
                               `1` = "Democratic Women")) %>%
  dplyr::rename(estimate = AME,
                conf.low = lower,
                conf.high = upper,
                term = factor) %>%
  dplyr::select(-party, -SE, -z, -p) %>%
  mutate(term = str_sub(term, 10)) -> data

dwplot(data) +
  theme_bw() +
  steph_theme +
  xlab("") +
  theme(legend.justification=c(.02, .993), legend.position=c(.02, .99),
        legend.title = element_blank(), legend.background =
          element_rect(color="gray90")) +
  scale_colour_grey(start = .3, end = .7,
                    name = "model") +
  geom_vline(xintercept = 0, colour = "grey60", linetype = 2) -> Fig6

Fig6

# ---------------------------------------------------------------------------------------------------------------
# Tables 1 & 2
# ---------------------------------------------------------------------------------------------------------------

# STUDY 1:
df <- readRDS("Qualtrics_2019.RDS")
df$CLINTON <- car::recode(df$TREATMENT, "'Alexa'= 0; 'Trump'=0; 'MeToo'=0; 'Marches'=0; 'Clinton'= 1")
df$TRUMP <- car::recode(df$TREATMENT, "'Alexa'= 0; 'Trump'=1; 'MeToo'=0; 'Marches'=0; 'Clinton'= 0")
df$METOO <- car::recode(df$TREATMENT, "'Alexa'= 0; 'Trump'=0; 'MeToo'=1; 'Marches'=0; 'Clinton'= 0")
df$MARCH <- car::recode(df$TREATMENT, "'Alexa'= 0; 'Trump'=0; 'MeToo'=0; 'Marches'=1; 'Clinton'= 0")
df$Engage <- df$futureengage
df$anxiety <- df$fear

mix_mod <- '
Engage ~ b1 * anxiety + b2 * anger + b3 * hope + c1 * TRUMP + c2 * CLINTON + c3 * MARCH + c4 * METOO
anxiety ~ a1 * TRUMP + a99 * CLINTON + a98 * MARCH + a97 * METOO
anger ~ a2 * TRUMP + a299 * CLINTON + a298 * MARCH + a297 * METOO
hope ~ a3 * TRUMP + a399 * CLINTON + a398 * MARCH + a397 * METOO
indirTRUMPanxiety := a1 * b1
indirTRUMPanger := a2 * b2
indirTRUMPhope := a3 * b3
indirCLINTONanxiety := a99 * b1
indirCLINTONanger := a299 * b2
indirCLINTONhope := a399 * b3
'
summary(lavaan:::sem(mix_mod, data=df, std.lv=TRUE))

# STUDY 2:
df_TESS <- readRDS("NORC.RDS")
df_TESS$CLINTON <- car::recode(df_TESS$TREATMENT, "'Alexa'= 0; 'Trump'=0; 'MeToo'=0; 'March'=0; 'Clinton'= 1")
df_TESS$TRUMP <- car::recode(df_TESS$TREATMENT, "'Alexa'= 0; 'Trump'=1; 'MeToo'=0; 'March'=0; 'Clinton'= 0")
df_TESS$METOO <- car::recode(df_TESS$TREATMENT, "'Alexa'= 0; 'Trump'=0; 'MeToo'=1; 'March'=0; 'Clinton'= 0")
df_TESS$MARCH <- car::recode(df_TESS$TREATMENT, "'Alexa'= 0; 'Trump'=0; 'MeToo'=0; 'March'=1; 'Clinton'= 0")
df_TESS$Engage <- df_TESS$futureengage
df_TESS$anxiety <- df_TESS$fear

mix_mod <- '
Engage ~ b1 * anxiety + b2 * anger + b3 * hope + c1 * TRUMP + c2 * CLINTON + c3 * MARCH + c4 * METOO
anxiety ~ a1 * TRUMP + a99 * CLINTON + a98 * MARCH + a97 * METOO
anger ~ a2 * TRUMP + a299 * CLINTON + a298 * MARCH + a297 * METOO
hope ~ a3 * TRUMP + a399 * CLINTON + a398 * MARCH + a397 * METOO
indirTRUMPanxiety := a1 * b1
indirTRUMPanger := a2 * b2
indirTRUMPhope := a3 * b3
indirCLINTONanxiety := a99 * b1
indirCLINTONanger := a299 * b2
indirCLINTONhope := a399 * b3
'

summary(lavaan:::sem(mix_mod, data=df_TESS, std.lv=TRUE))

# ---------------------------------------------------------------------------------------------------------------
# Table 3: Treatment effects among compliers, CCES 2020
# ---------------------------------------------------------------------------------------------------------------
df_cces <- readRDS("CES2020.RDS")

# NonCompliers Subset
df_cces %>%
  filter(complier == 1) -> df_compliers
mod2 <- lm(futureengage_novote ~ TREATMENT, data = df_compliers)
export_summs(mod2)

# ---------------------------------------------------------------------------------------------------------------
# In Text Material: Compliance and Engagement

mod <- lm(complier ~ Gender_Party, data = df_cces)
summary(mod)
emm <- emmeans(mod, ~ Gender_Party) 
summary(contrast(emm, method = "revpairwise", 
                 interaction = TRUE, adjust = "none"))
emm

# By Treatment:
mod <- lm(complier ~ Gender_Party * TREATMENT, data = df_cces)
emm <- emmeans(mod, ~ Gender_Party | TREATMENT) 
summary(contrast(emm, method = "revpairwise", 
                 interaction = TRUE, adjust = "none"))
emm

# Did one group write more?

mod <- lm(word_count ~ Gender_Party, data = df_cces)
summary(mod)
emm <- emmeans(mod, ~ Gender_Party) 
summary(contrast(emm, method = "revpairwise", 
                 interaction = TRUE, adjust = "none"))
emm

# By Treatment

mod <- lm(word_count ~ Gender_Party * TREATMENT, data = df_cces)
summary(mod)
emm <- emmeans(mod, ~ Gender_Party | TREATMENT) 
summary(contrast(emm, method = "revpairwise", 
                 interaction = TRUE, adjust = "none"))
emm

